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ABSTRACT 


This  paper  presents  a  new  algorithm  to  solve  a  network  problem 
with  equal  flow  side  constraints.  The  proposed  solution  technique  is 
motivated  by  the  desire  to  exploit  the  special  structure  of  the  side 
constraints  and  to  maintain  as  much  of  the  characteristics  of  pure 
network  problems  as  possible.^  Not  only  has  specialized  software  for  the 
efficient  solution  of  ptire  networks  been  developed,  but  the  same  compu¬ 
tational  efficacies  lend  themselves  to  the  solution  of  sequences  of 
minimum  cost  network  flow  problems  by  using  reoptimization  procedures. 

Our  solution  technique  for  the  equal  flow  problem  consists  of  solving 
two  sequences  of  pure  network  problems.  One  sequence  yields  tighter 
lower  bounds  on  the  optimal  value  by  considering  the  Lagrangean  relax¬ 
ation  of  the  equal  flow  problem  in  which  the  side  constraints  are  dualized 
The  second  sequence  yields  upper  bounds  on  the  optimal  value  for  the 
problem  and  maintains  a  feasible  solution  at  all  times.  This  sequence  is 
obtained  by  considering  a  reformulation  of  the  equal  flow  problem  based 
on  parametric  changes  in  the  requirements  vector.  The  procedure  has  the 
added  attractive  feature  that  it  provides  a  feasible  solution  which  is 
known  to  be  within  a  percentage  of  the  optimal  at  all  times.  As  such, 
the  algorithm  terminates  when  a  solution  with  a  prespecified  tolerance 
on  the  objective  function  value  is  obtained.  On  NETGEN  problems, 
using  the  first  150  arcs  to  form  75  equal  flow  side  constraints,  we  found 
that  the  new  algorithm  is  approximately  3  times  faster  than  existing 
techniques  and  requires  only  50%  of  the  storage. 
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I.  INTRODUCTION 


This  paper  presents  a  new  technique  to  solve  the  equal  flow  problem. 
This  problem  is  easily  conceptualized  as  a  minimal  cost  network  flow 
problem  with  additional  constraints  on  certain  pairs  of  arcs.  Specifically, 
given  pairs  of  arcs  are  required  to  take  on  the  same  value.  Applications 
of  this  model  include  crew  scheduling  [6],  estimating  driver  costs  for 
transit  operations  [28],  and  the  two  duty  period  scheduling  problem  [25]. 

The  equal  flow  problem  may  be  solved  using  a  specialization  of  the  simplex 
method  for  networks  with  side  constraints.  However,  by  exploiting  the 
special  structure  of  the  side  constraints,  we  have  developed  a  new  algo¬ 
rithm  which  results  in  a  decrease  in  both  computer  storage  and  compu¬ 
tation  time. 

It  is  well  documented  that  pure  network  problems  can  be  solved  from 
fifty  to  one  hundred  times  faster  using  specialized  primal  simplex  soft¬ 
ware  as  compared  to  general  linear  programming  systems.  Motivated  by  this 
great  advantage,  our  procedure  solves  the  equal  flow  problem  as  a  sequence 
of  pure  network  problems  and  totally  eliminates  the  need  to  deal  with  a 
basis  matrix. 

1.1  Problem  Description 

The  equal  flow  problem  is  defined  on  a  network  represented  by  an  (m,n) 
node-arc  incidence  matrix,  A,  in  which  K  pairs  of  arcs  are  identified  and 
required  to  have  equal  flow.  Mathematically,  this  is  expressed  as: 


i  r  i  *  iJ  JLid 


Minimize  cx 


s.t.  Ax  ■  b 

*k  "  *k+K  *  k  «  1,  ....  K 

0  1  x  i  u 

where,  c  is  a  1  x  n  vector  of  unit  costs,  b  is  an  m  x  1  vector  of  node 
requirements,  0  is  an  n  x  1  vector  of  zeroes,  x  is  an  n  x  1  vector  of 
decision  variables,  and  u  is  an  n  x  1  vector  of  upper  bounds.  The  above 
definition,  henceforth  referred  to  as  Pi,  assumes  that  the  first  2K  arcs 
appear  in  the  equal  flow  constraints.  This  assumption  is  in  no  way 
restrictive  since,  by  rearranging  the  order  of  the  arcs,  any  equal  flow 
problem  with  K  pairs  can  be  expressed  in  the  above  form.  Note  that  the  K 
pairs  of  arcs  are  mutually  exclusive,  i.e.  an  arc  appears  in  at  most 
one  side  constraint. 

1.2  Survey  of  Related  Literature 

In  1961,  Charnes  and  Cooper  [7]  presented  a  specialized  algorithm 
for  the  model: 

Minimize  cx 
s.t.  Ax  ■  b 

Cx  -  d 
x  >  0, 

where  A  and  C  are  some  general  matrices  but  A  has  some  favored  structure. 
Their  algorithm,  called  the  double  reverse  method,  takes  advantage  of 
the  special  structure  of  the  matrix  A.  Variations  of  this  algorithm  may 
be  found  in  [2,  8,  10,  L5,  19,  24] .  Specializations  for  multicommodity 


problems  may  be  found  in  [14,  20,  21]. 

In  1980,  Shepardson  ,and  Marscen  [25]  showed  chat  the  two  duty  period 
scheduling  problem  can  be  reformulated  as  a  single  duty  period  scheduling 
problem  with  equal  flow  side  constraints.  They  obtain  a  Lagrangean  dual 
for  this  equal  flow  problem,  by  dualizing  with  respect  to  the  equal  flow 
side  constraints.  This  Lagrangean  dual  is  maximized  using  the  subgradient 
optimization  technique.  In  1984,  Turnquist  and  Malandraki  [28]  modeled  the 
problem  of  estimating  driver  costs  for  transit  operations  as  an  integer 
equal  flow  problem.  They  obtain  a  Lagrangean  dual  for  their  problem,  by 
dualizing  with  respect  to  the  side  constraints.  Their  algorithm  is  a 
slight  modification  of  the  subgradient  optimization  technique.  They 
perform  a  line  search  between  two  successive  solutions  obtained  during 
the  subgradient  optimization  process. 

Beck,  Lasdon,  and  Engquist  [5]  transformed  the  equal  flow  problem 
into  a  quadratic  programming  problem  which  has  a  penalty  for  violating 
the  equal  flow  constraints.  They  solved  this  nonlinear  programming 
problem  ueing  the  Fletcher-Reeves  conjugate  gradient  method  [9] ,  a 
successive  linear  programming  code  [13 ] ,  and  a  convex  simplex  code. 

If  the  penalty  is  sufficiently  large,  this  approach  is  guaranteed  to 
converge  to  the  optimal  solution  of  the  equal  flow  problem. 

1.3  Objective  of  the  Investigation 

The  objective  of  this  investigation  is  to  develop  and  computationally 
test  a  new  algorithm  for  the  equal  flow  problem.  This  algorithm  utilizes 
the  subgradient  optimization  technique  and  is  based  on  the  relaxation/ 
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restriction  procedure  proposed  by  Glover,  Glover,  and  Martinson  [11]  for 
a  generalized  network  model  with  special  side  constraints.  We  establish 
that  the  equal  flow  problem  may  be  solved  as  two  sequences  of  pure 
network  problems,  one  sequence  corresponds  to  computing  a  lower  bound 
while  the  other  corresponds  to  computing  an  upper  bound.  In  the  limit, 
both  bounds  will  converge  to  the  optimal  objective  value.  Our  implementation 
terminates  when  the  difference  between  the  bounds  is  within  a  prespecified 
tolerance. 

The  subgradient  optimization  technique  requires  the  computation  of 
subgradients,  choice  of  appropriate  step  sizes,  and  the  application  of 
a  projection  operation.  We  show  that  the  subgradients  for  the  upper 
bound  can  be  computed  using  the  optimal  dual  variables  obtained  by 
solving  pure  network  problems.  We  also  develop  theoretical  results 
that  yield  an  easy  implementation  of  the  projection  operation.  The 
step  sizes  selected  are  a  modification  of  the  ones  proposed  by  Polyak 
[23].  For  this  choice  of  step  sizes,  we  prove  that  our  algorithm  must 
necessarily  obtain  an  iterate  at  which  the  objective  value  is  arbitrari¬ 
ly  close  to  the  optimal  objective  value.  In  a  computational  study, 
comparing  our  code  with  a  code  chat  is  designed  to  solve  network  problems 
with  side  constraints,  we  found  that  the  new  code  runs  approximately  3 
times  faster  and  requires  50Z  less  core  storage. 
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II.  THE  SUBGRADIENT  ALGORITHM 


The  Subgradient  Algorithm  was  first  introduced  by  Shor  [27]  and  is 
a  general  procedure  for  solving  nonlinear  programming  problems.  It  may 
be  viewed  as  a  generalization  of  the  steepest  descent  (ascent)  method 
for  convex  (concave)  problems  in  which  the  gradient  may  not  exist 
everywhere.  The  subgradient  is  simply  substituted  in  place  of  the 
gradient  for  those  points  for  which  the  gradient  does  not  exist.  When 
this  occurs,  the  algorithm  may  move  to  a  point  with  objective  value 
worse  than  the  current  point.  Hence,  the  objective  function  does  not 
necessarily  improve  at  each  iteration  and  consequently  the  convergence 
results  of  Zangwill  [29]  do  not  apply.  Remarkably  though,  under  fairly 
minor  conditions  on  the  step  size,  convergence  can  be  guaranteed. 

Let  the  nonlinear  program  PO  be  given  by: 

Minimize  f(y) 
s.t.  y  z  G 

where  f  is  a  real  valued  function  that  is  convex  over  the  compact,  con¬ 
vex,  and  nonempty  set  G.  A  vector  n  will  be  called  a  sub gradient  of  f 
at  y  if  f(y)  -  f(y)  >  n(y  -  y)  for  all  y  e  G.  For  any  y  e  G,  we  denote 
the  set  of  all  subgradients  of  f  at  y  by  3f(y).  The  subgradient  algorithm 
makes  use  of  an  operation  called  the  projection  operation.  The  projection 
of  a  point  x  onto  G,  denoted  by  P[x],  is  defined  to  be  the  unique  point  y  e  G 
that  is  nearest  to  x  with  respect  to  the  Euclidean  norm.  Using  the  projec¬ 
tion  operation,  we  now  present  the  subgradient  algorithm  in  its  most 
general  form. 
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A LG  1  SUBGRAVIENT  OPTIMIZATION  ALGORITHM 


Step  0  (Initialization) 

Let  y^  be  any  element  of  G,  select  a  set  of  step  s-izes 
Sgfs^,S2,  •  •• ,  and  set  i  0. 

Step  1  (Find  Subgradient) 

Let  n±  e  df(yp.  If  *  0,  then  terminate  with  y^  optimal. 

Step  2  (Move  to  New  Point) 

Set  y^  P[yi  -  s^] ,  set  i  i  +  1  and  return  to  step  1. 

Various  proposals  have  been  offered  for  the  selection  of  the  step 
sizes.  Three  general  schema  which  have  been  suggested  are: 


ii) 


iii) 


Xi(f(yi)-f*) 

IlnJI2 


o  <  \  <  2, 


where  f*  is  the  optimal  value  of  f  over  G.  If  the  constants,  A^'s, 
satisfy  the  following  conditions: 

A.  >  0,  all  i;  lim  A.  *0;  and  £  A.  »  », 

X  ±*  as  1  1 

then  the  convergence  of  the  algorithm  is  guaranteed  using  (i)  or  (ii) 
(see  Goff in  [12],  Helgason  [17],  Kennington  and  Helgason  [21]).  For  the 
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upper  bounds,  we  use  a  modification  of  the  third  step  size.  The 
following  result  is  available  for  this  scheme. 

Proposition  1  (Polyak  [23]) 

Let  f  be  a  real  valued  convex  function  over  the  compact,  convex, 
and  nonempty  set  G.  Also,  let  f*  be  the  minimum  of  f  and  ||n^||  <  C  for 
all  i  and  some  constant  C.  Then  there  exists  a  y*  s  G  with  f(y*)  -►  f*, 
if  scheme  (iii)  is  used. 

Note  that  in  (iii),  f*  is  the  optimal  value  of  f  over  G.  Since  the 
optimal  objective  is  unknown  before  solving  the  problem,  we  use  a  lower 
bound  on  f*  in  our  implementation. 


III.  THE  LOWER  BOUND 


Recall  thac  Che  equal  flow  problem,  which  vt  denote  by  PI,  is 
given  by: 

Minimize  cx 
s.c.  Ax  -  b 

*k  "  *K+k  *  k  "  1,***’K 
0  <  x  <  u. 

In  our  algorithm  for  PI,  lower  bounds  on  the  optimal  objective  of  PI  are 
used  for  step  sizes  and  for  termination.  In  this  section,  we  describe  a 
procedure  to  obtain  these  lower  bounds. 

Consider  the  following  Lagrangean  dual  for  PI,  which  we  shall 
refer  to  as  Dl:  _ 

Maximize  h(w)  ,  where  w  *  [w  , ...,w  ]  z  R\  and 
K 

h(w)  ■  Min{cx  +  £  :  ^  ^  ®  1  x  1 

k-1  ' 

Proposition  2  (Shetty  [26]) 

Let  x  he  a  feasible  solution  to  PI  and  let  w  *  [w^,...,w^]  be 
a  feasible  solution  to  Dl.  Then  cx  >  h(w) . 

Proposition  3  (Bazaraa  and  Shetty  [4]) 

If  PI  has  a  minimum,  then  the  optimal  objectives  for  PI  and  Dl 
are  equal. 

As  a  consequence  of  Propositions  2  and  3,  we  may  solve  Dl  to 
obtain  a  lower  bound.  We  will  now  show  chat  Dl  may  be  solved  using  the 
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VI .  COMPUTATIONAL  EXPERIMENTATION 


This  section  describes  the  computer  implementation,  EQFLO,  and 
testing  of  our  algorithm  for  the  equal  flow  problem.  The  algorithm 
was  tested  on  a  set  of  35  test  problems  randomly  generated  using 
NETGEN  [22].  Computation  times  are  compared  with  those  of  NETSIDE  [2], 
a  general  purpose  code  for  network  problems  with  side  constraints. 

Both  NETSIDE  and  EQFLO  are  written  in  standard  FORTRAN  for  an  incore 
implementation  and  have  not  been  tailored  to  either  the  machine  or 
FORTRAN  compiler  used  for  testing. 

6.1  Description  of  the  Computer  Codes 

NETSIDE  was  developed  by  Barr,  Farhangian,  and  Kennington  at 
Southern  Methodist  University,  Dallas,  Texas.  Designed  to  solve  network 
problems  with  side  constraints,  it  used  a  specialization  of  the  revised 
simplex  method  known  as  the  primal  partitioning  algorithm  [15],  The 
basis  inverse  is  maintained  as  a  rooted  spanning  tree  and  a  working  basis 
inverse  in  product  form.  The  reinversion  routine  is  a  modification  of  the 
work  of  Hellerman  and  Rarick  [18]  and  uses  the  "spike  swapping  theory" 
of  Helgason  and  Kennington  [16].  The  initial  working  basis  consists  of  a 
combination  of  artificial  and  slack  variables.  The  working  basis  is 
reinverted  every  50  iterations.  The  pricing  routine  uses  a  candidate  list 
of  size  10  with  a  block  size  of  400.  Both  pricing  and  pivot  tolerance  are 
l.E-6. 

EQFLO  is  our  iaplementaiton  of  ALG  4,  and  makes  use  of  MODFLO  [1]  to 
solve  pure  network  subproblems.  MODFLO  is  a  set  of  subroutines  which  may 
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Seep  2  (Compute  Upper  Bounds) 

2a.  Set  T  +  0,  set  IFLAG  *■  1. 

2b.  Call  Alg  3  (steps  2  and  3a) 

2c.  Set  R  R+l. 

If  R  <  ITERU,  then  go  to  2b. 

Step  3 

Set  R  «-  0. 

Set  p  -  pQ. 

Go  to  1. 

In  the  above  algorithm,  IFLAG  is  used  in  obtaining  a  starting  y  from  the 
solutions  in  la.  The  bases  used  in  steps  1  and  2  are  generated  from  the 
optimal  bases  obtained  in  the  previous  iterations. 
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V.  THE  ALGORITHM 


In  this  section,  we  present  our  new  algorithm  for  solving  the  equal 
flow  problem.  Let  ITERL  denote  the  number  of  iterations  spent  in  step  la 
in  computing  the  lower  bound  before  returning  to  the  upper  bound,  and 
ITERU  denote  the  number  of  iterations  spent  in  computing  the  upper  bound 
before  returning  to  the  lower  bound.  Also,  let  T  denote  the  iteration 
count  for  the  lower  bound,  R  denote  the  iteration  count  for  upper  bound, 
and  Pq  denote  the  initial  step  size  for  the  lower  bound. 

ALG  4  SUBGRAVIENT  OPTIMIZATION  ALGORITHM 
FOR  THE  EQUAL  FLOW  PROBLEM 

Step  0  (Initialization) 

Initialize  ITERL,  ITERU,  REREQ,  XQ,  pQ,  and  tolerance  s. 

Sec  T  0,  sec  Q  •*-  0,  set  R  0,  sec  w  0,  and  sec  IFLAG  *■  0. 

Sec  UBND  -*-  +  •,  set  LBND  ■*--•,  set  p  pQ,  and  set  *■  for  k»l,.... 
Set  u  (minCu^Uj^) , . . .  ,min(u^,  u2R) ) . 

Step  1  (Compute  Lower  Bounds) 

la.  Call  ALG  2  (steps  2  and  3a) . 

lb.  Sec  T  -  T+l 

If  T  •**  ITERL,  then  go  to  1. 

lc.  (Initialize  y) 

If  IFLAG  A  0,  'then  go  to  2;  otherwise, 

sec  y  *■  (minCu^,  (x^+^.+i)/2) , . . .  ,min(u^,,  2K)/2)]* 
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6  2  6 

-  - - —  +  g*( - a) 

(2Ai  -  Xp  2-6  2-6 

2  6 

<  M6  +  g*(  -  — - -  a),  where  M  is  a  constant  less  than 

^  “  2-6  2-6 

(2e-e2) 

This  completes  the  proof  of  Proposition  12. 
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We  can  choose  an  inceger  N  large  enough  chac 


(g*  ~  i)6 


Adding  together  the  inequalities  obtained  from  (1)  by  letting  i  take 


on  all  values  from  1  to  N,  we  obtain 


2  -  N(g*  -  g)6 

II  yN  -  y*  II  <  II  y±  -  y*  II  -  - 2 -  <  °* 

c 


a  contradiction.  This  justifies  our  assertion. 


By  simplifying  our  assertion  further  we  get, 


6  g  A.  2g* 

g(y  )  <  - - - +  - 

(2Ai  -  X1)  (2  -  Xp  (2  -  A.) 


a  g*A  2g* 

- = —  + - 

(2  -  Ap  (2  -  Ai) 


- - 7—  +  g*(l  + - - -  (1  -  a)) 

(2Ai  -  Ap  (2  -  Xi) 


5  3 

- -  +  g*(i  + -  (1  -  a)) 

(2A  -  Ap  2-8 
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Suppose  chat  for  all  i. 


Ai 1 

g<yt)  > - 5 - 


2A.  g* 


(X  -  2A  )  Ck.  -  2X  . )  a  -  2A . ) 

i  1  ii  i  i' 


,  where 


6  >  0  is  given.  Let  y*  e  S  be  an  optimal  point.  By  Proposition  11, 

^?(8(y,)-I)2  2Ai(g(y i)”i)  ni(y*-y1) 


yi+ry*U2i  II  yi-y*ll2  + 


i  ”  w  i 


II  n. 


II 


Since  e  3g(yi),  ni<y*  -  y± )  <  g*  -  g(yi). 
Thus, 


2  _  „  ..  _„2  1  x 


yi+i-y*ll  <  II  y±- y*ll  + 


^Cg(y,)-I)2  2Xi(g(yi)-i)(g*  -  g^)) 


II  n± 


y±  -  y*  ||2  +  (g(yi)-i) 


g(yt) (Xt 


ni 


2X^)-X^  g  +  2Ai  g* 


Hill 


i  II  y±  -  y*  II2  - 


(sCy^  -  1)5 

ii  \  ii2 


i  II  y,  -  y 


-  -.1.2  <»*-!>* 


(1) 
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v- ■>  . 

.►  .•  J  .A.’  »-•  »■-  ■  «-'■  «-•  «- 


In  our  implementation,  we  use  g,  a  lower  bound  on  g,  in  place  of  g*. 

The  following  propositions  deoonscrate  chat  for  g's  close  to  g*,  our 
procedure  oust  necessarily  obtain  an  iterate  at  which  g  is  arbitrarily 
close  to  g*. 

Proposition  11  (Kennington  and  Helgason  [21]) 

If  *  0, 

II  yi+1  -  y±  II  2  <  II  y±  -  y||2+  s2  II  nil|2+  Zs^Cy  -  y±)  for  any  y  e  S 
and  step  size  s^. 

Proposition  12 

Let  g*  be  the  optimal  value,  of  g ,  and  also  let 
i)  xx  g*  <  g  <  g*,  0  <  a  <  1, 
ii)  si  -  XjL(g(yi)  -  i)/||  |1  2  ,  and 

iii)  0<e<A^<8<2  for  all  i. 

If  there  is  a  constant  C  such  that  ||  rij.ll  <  C  for  all  i,  then  there 

2  8 

exists  some  i  such  that  gCy^  <  M6  +  g*(  2-3"  “  “Jig”  ®)  for  ax*y  >  0 
and  for  some  constant  M. 

Proof 

First,  we  assert  that  there  is  some  i  such  that 

A \  i  2At  g*  5 

g(y.)  <  - 3 - 5 - 5 -  ,  for  any  5  >  0. 

1  (I’  -  l\J  (A*  -  2A.)  (X‘  -  2At) 
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Set  UBND  cx 
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where  L,M  are  integers  and  0  <  L,M  <  K.  Then 

P(y)  -  y*  -  (0,....0,yL+1,...,?L+tl,GL+iH.1 . uR)  is 

a  projection  of  y  on  S. 

Following  a  description  of  the  terminology  used,  our  algorithm  for 
obtaining  an  upper  bound  for  PI  is  presented  below.  Let  RFREQ  denote  the 
frequency  at  which  the  constant  X  in  step  size  (iii)  is  reset  to  its  initial 
value,  X  ,  LBND  denote  a  lower  bound  on  PI,  UBND  denote  an  upper  bound  on 
PI,  P  denote  the  projection  routine  described  in  Proposition  10,  e  denote 
the  termination  tolerance  and  Q  denote  the  iteration  count  for  the 
upper  bound. 


A LG  3  UPPER  BOUND  ALGORITHM 

Step  1  (Initialization) 

Choose  y  e  S. 

Initialize  LBND,  RFREQ,  e,  and  XQ. 

Set  u  -*•  (min(u1(uK+1) ,...  .minCu^.u^)) . 

Set  Q  0,  sec  ^  «■  X^  for  k»l,...,K. 

Step  2  (Find  Subgradient  and  Step  Size) 

w  * 

For  allocation  y,  let  x  and  it,  respectively,  be  the  vectors  of 

optimal  primal  and  dual  variables  for 

Min{cx:  Ax  *  o,  0  <  x  <  u}.  Construct  x  from  x,y. 
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for  P4.  Suppose  the  arc  corresponding  to  j  has  "Prom"  node  j ^  and  "To" 
node  jy  That  is,  arc  j  is  the  ordered  pair  (j1>j2)*  Then  we  define 
FROM(j)  -  ^  and  TO(j)  *  j2*  Using  this  notation,  the  following  propo¬ 
sition  gives  the  required  formulae: 


Proposition  8  (Shetty  [26]) 

Let  ir  be  the  vector  of  optimal  dual  variables  for  P4.  Then 
[*>  v,  ,vv^,...,vv,v,v]  with  -  -PFR0M(j)  +  *T0(j)+  cj’  j"1*---’2*’ 


'1* ’K+l’ "  ”  *  K’  2K 
are  optimal  duals  for  ?3. 


FROM(j) 


We  now  present  two  propositions  that  justify  the  projection  routine 
used  for  the  upper  bound.  The  proofs  may  be  found  in  Kennington  and 
Helgason  [21]  and  Shetty  [26]. 


Proposition  9 

Let  S  be  a  nonempty,  convex  set  and  y  t  S.  Then  y*  e  S  is  a 
projection  of  y  on  to  S  if  (y  -  y*)  (y  -  y*)  <  0  for  all  y  e  S. 


Proposition  10 

Uty*  ^1*^2 . £  rK  with 

yk  <  0  for  k  ■  1, . . .  ,L 

0  <  y^  <  uk,  -  min^,  uR+k)  for  k  •  L+1,...,L+M 
yk  >  Uj^  for  k  ■  L+Mfl, . . .  ,K 
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Minimize  cx 


s.t.  Ax  *  b 

(IT) 

X1  '  yl 

<V 

*K+1  "  yl 

• 

(W 

e 

• 

“  yK 

(V 

X2K  "  y2K 

<v2k> 

0  <  x  <  u 

(y) . 

"  ■  (vi  +  ,W---vk  + 

v2K>  is  a  subgradient jjf  g  at  y  *  (y^  • . .  ,yK>  . 

As  a  result  of  Proposition  7,  a  subgradient  at  any  given  point 


y  ■  (ylf...,yK)  £  S  required  in  our  specialization  of  ALG  1  can  be 


obtained  by  solving  Min{cx:  Ax  *  b,  x^  *  "  yK’  0  -  x  -  u^»  w^ich 


we  shell  refer  to  as  P3.  After  substituting  x.  *  y.,...,x,v  ■  y  ,  in 

X  X  Xiv  iv 


Ax-*  b,  we  obtain  a  pure  network  problem,  which  we  shell  refer  to  as 
P4  and  is  given  below: 

Minimize  cx 
s  •  c  • 


AA  A 

Ax  ■  b 


-  ,  A  A 

0  <  x  <  u. 

To  apply  ALG  1,  we  need  a  procedure  for  constructing  the  optimal  dual 
variables  (vi»vK+1» • • • »vk»v2k^  for  ?3  froa  ch*  op*1®*1  dual  variables 


IV.  THE  UPPER  BOUND 


An  alternate  formulation  of  PI,  which  will  be  referred  to  as  P2, 
is  as  follows: 

Minimize  g(y) 
s.t.  y  £  S 

where  for  any  vector  y  *  [y^,...,yK], 

g(y)  -  Min{cx:  Ax  -  b,  0  <  x  <  u,  -  ^  for  a11 


S  -  {y:  0  <  yk  <  min(uk,uK+k)  for  all  It}. 

Clearly,  PI  and  P2  are  equivalent.  That  is,  given  an  optimum  for  one, 
we  can  construes  an  optimum  for  the  other.  We  will  now  show  that  P2 
is  a  special  case  of  the  nonlinear  program  PO  and  may  be  solved  using 
the  Subgradient  Optimization  Algorithm,  ALG  1. 

Proposition  6  (Shetty  [261) 

The  real  valued  function  g  is  piece-wise  linear  convex  over  the 
compact,  convex  and  nonempty  sec  S. 

To  apply  the  subgradienc  algorithm,  we  need  a  procedure  for  obtaining 
a  subgradient  of  g  at  a  point  y.  The  following  proposition  shows  that 
the  dual  variables  may  be  used  to  construct  a  subgradient. 


Proposition  7  (Shettv  f26]) 

Let  <1T»v1»vK+1»‘**’vK>v2K'y) 


be  the  optimal  dual  variables  for 
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Step  2  (Find  Subgradient) 

Let  x  -  [x, ]  solve 
1  n 

K 

h(w)  «  Min{cx  +  l  -  x^) :  Ax  -  b,  0  <  x  <  u}. 

k"l 

Set  LBND  -  h(w) . 


If  (UBND  -  LBND)  <  e(UBND)  then  terminate; 
otherwise,  set  d  [(x^  -  » •  •  • »  (*K  ”  . 

Step  3  (Move  to  a  New  Point) 

3a.  Set  w  w  +  pd,  set  p  **•  p/2. 

3b.  Go  to  2. 


m 


v 


■ 


- .  *1 


I 


subgradient  optimization  technique  for  concave  functions.  This  tech¬ 
nique  is  similar  to  ALC  1  with  a  modification.  Let  Po»P1»P2»“* 
denote  a  sequence  of  step  sizes  and  let  di  t  ShCw^ .  Then  step  2  is 
replaced  by: 

Step  2  (Move  to  New  Point) 

Set  wi+1  w^  +  p^d^,  set  i  i+1  and  return  to  step  1. 

To  use  this  algorithm  h(w)  must  be  concave,  and  we  need  a  means 
of  generating  subgradients.  These  two  results  follow: 


B 


Proposition  4  (Shetty  f 26 1 ) 

ir 

The  real  valued  function  h  is  concave  over  R  , 

Proposition  5  (Shettv  f26P 

For  a  given  w,  let  x  be  an  optimal  solution  to 
K 

Min{cx  +  £  w^(x^  ~  :  ®  1  x  1  u^' 

k«l 

Then  d  ■  [(x^  -  x^,^) , . . . ,  (x^  -  x2K>  1  is  a  subgradient  of  h  at  w. 

We  used  scheme  (i)  for  step  sizes.  Let  UBND  denote  an  upper  bound  and 
assume  that  the  optimal  objective  value  is  positive.  Our  algorithm  for 
obtaining  lower  bounds  is  presented  below: 

A LG  2  LOWER  SOUNV  ALGORITHM 

Step  1  (Initialization) 

Initialize  UBHD,  step  size  p,  and  tolerance  e. 

Sec  w  *■  0. 
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be  used  to  solve  a  network  problem  as  well  as  reoptimize  after  problem 
data  changes.  Based  on  NETFLO  [21],  this  code  allows  the  user  to  change 
costs,  bounds  and/or  requirements  for  a  network  problem  and  reoptimize. 

The  tuning  parameters  used  in  all  runs  were  as  follows:  ITERL  *  5, 

ITERU  *  10,  REFREQ  -  5,  XQ  -  0.75,  pQ  -  0.01,  and  £  =»  0.1.  MODFLO  [1] 
is  used  to  reoptimize  after  each  change  to  either  the  costs  or  right- 
hand-sides. 

6.2  The  Test  Problems 

The  program  NETGEN,  a  generator  for  large-scale  network  test 
problems,  was  used  to  generate  35  test  problems.  The  parameters  used 
to  generate  these  problems  are  described  in  Klingman,  Napier,  and  Stutz  [22]. 
The  test  problems  have  between  200  and  1500  nodes,  and  1500  and  7000  arcs. 

For  each  problem,  the  first  150  arcs  were  paired  to  form  equal  flow  sides 
constraints.  The  characteristics  of  these  test  problems  are  listed  in 
Table  1. 

Our  algorithm  requires  upper  bounds  on  all  equal  flow  arcs.  Though 
NETGEN  generates  bounds  on  some  of  these  arcs,  there  were  others  with  no 
upper  bounds.  We  3et  the  maximum  of  all  supplies  and  demands  to  be  the 
upper  bounds  on  such  arcs.  These  bounds  were  acceptable  since  the  optimal 
solutions  obtained  for  our  pure  network  problems  were  the  same  as  the  ones 
listed  in  NETGEN.  Furthermore,  for  all  35  test  problems  the  first  150  arcs 
were  _sed  to  form  75  pairs  of  equal  flow  side  constraints.  We  were  unable 
to  experiment  with  more  than  75  pairs  due  to  a  core  storage  limitation  of 
301K.  NETSIDE  required  approximately  300K  octal  words  of  storage  for 
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problems  28  through  35  with  75  pairs  and  any  further  increase  in  the 
number  of  pairs  would  exceed  the  storage  limitation. 


set  Q  Q+l. 

Note  that  this  modification  results  in  different  step  sizes  for  different 
equal  flow  pairs.  The  details  of  all  runs  are  given  in  Tables  2  and  3. 
The  times  are  in  CPU  seconds  and  exclude  input  and  output. 

The  value  0.01  used  for  p^  worked  well  for  all  test  problems  except 
problem  number  11.  This  problem  experienced  difficulties  in  converging 
within  10%  of  the  optimal.  However,  the  problem  did  converge  within  10% 
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of  Che  optimum  when  we  changed  Pq  to  values  between  3  and  10. 

The  computational  results  presented  in  Tables  2  and  3  are  summarized 
in  Table  4.  Letting  T(ALG)  denote  the  CPU  time  required  to  solve  the  35 
test  problems  using  code  ALG,  the  relationship  is  given  below: 

T( NETS IDE)  -  2.54  (EQFLO) ,  same  step  size, 

T(NETSIDE)  *  3.00  (EQFLO),  different  step  sizes. 

Note  that  EQFLO  performs  better  as  the  problem  size  increases.  Although 
EQFLO  was  slightly  slower  than  NETSIDE  on  problems  1  to  10,  its  performance 
increased  substantially  on  problems  11  to  35.  In  particular,  EQFLO  ran 
approximately  5  to  6  times  faster  than  NETSIDE  on  problems  28  through  35 
and  these  problems  are  fairly  large.  We  expect  EQFLO  to  perform  even 
better  on  much  larger  problems.  This  is  attributable  to  the  fact  that 
the  time  for  pricing  and  updating  increases  dramatically  for  NETSIDE 
with  an  increase  in  the  size  of  the  network.,  whereas,  the  time  increase 
should  be  relatively  small  for  EQFLO  because  the  above  operations  are 
performed  very  efficiently  using  labelling  procedures  on  the  rooted 
spanning  tree. 

The  75  side  constraints  made  the  problems  approximately  three  times 
harder.  That  is,  the  pure  networks  were  solved  in  693  seconds  while  it 
required  1973  seconds  to  solve  the  equal  flow  problem.  Klingman, 

Napier,  and  Stutz  [22]  solved  the  same  35  pure  network  problems  in 
approximately  200  seconds  using  an  advance  start  on  a  CDC  6600  nr  rh» 
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University  of  Texas  at  Austin.  Richard  Barr's  best  tine  on  these 
35  test  problems  is  104  seconds  using  ARC  II  [3].  This  difference  in 
time  is  due  to  the  fact  that  EQFLO  is  a  real  code  (as  opposed  to  all¬ 
integer),  uses  an  all  artificial  start,  and  does  not  use  the  advanced 
data  structure  or  candidate  list  incorporated  in  ARC  II. 

These  35  problems  were  the  largest  that  could  be  solved  using 
NETSIDE  under  a  core  storage  limitation  of  301K  octal  words.  However, 
EQFLO  required  much  less  storage;  approximately  502  less  than  NETSIDE. 
This  additional  storage  for  NETSIDE  results  from  the  working  basis 
inverse  and  the  arrays  required  during  the  reinversion  process. 


VII.  SUMMARY  AND  CONCLUSIONS 


This  paper  presents  a  new  procedure  for  the  equal  flow  problem. 
Unlike  the  simplex  method  for  the  network  problem  with  side  constraints, 
this  new  procedure  does  not  require  a  working  basis.  We  have  showed 
that  using  the  subgradient  optimization  technique,  the  equal  flow 
problem  may  be  solved  as  two  sequences  of  pure  network  problems.  One 
sequence  corresponds  to  a  lower  bound  while  the  other  corresponds  to 
an  upper  bound.  In  the  lower  bound,  each  network  differs  from  the 
previous  one  in  that  the  cost  vector  has  changed.  In  the  upper  bound, 
each  network  differs  from  the  previous  one  in  that  the  right  hand  side 
has  changed.  While  solving  the  pure  network  problems  with  these 
changes  in  the  problem  data,  a  reoptimization  procedure  is  used  to 
obtain  a  good  starting  solution.  Our  technique  terminates  when  the 
difference  between  two  bounds  is  within  a  prespecified  tolerance. 

Subgradients  for  upper  bounds  are  computed  using  the  optimal 
dual  variables  obtained  by  solving  the  pure  network  problems.  The  sub¬ 
gradients  for  lower  bounds  are  the  difference  between  the  flows  on  the 
equal  flow  arcs,  obtained  while  solving  the  Lagrangean  relaxation.  The 
projection  operation  is  easily  implemented.  The  step  sizes  (i)  and 
(ill),  described  in  Section  II,  are  used  for  lower  and  upper  bounds, 
respectively.  For  these  step  sizes,  we  are  guaranteed  a  solution  at 
which  the  objective  value  is  arbitrarily  close  to  the  optimal  objective 
value. 

We  solved  all  test  problems  twice;  once  with  the  same  step  size  for 


all  equal  flow  pairs,  and  once  with  different  step  sizes  for  each  pair. 
The  tests  were  conducted  on  a  set  of  35  randomly  generated  problems  and 
a  comparison  was  made  with  NETSIDE,  a  code  that  is  designed  to  solve 
network  problems  with  side  constraints.  On  the  average,  our  code  ran 
approximately  3  times  faster.  Hovever,  it’s  performance  improved 
substantially  as  the  problem  size  increased.  The  new  algorithm  requires 
only  50%  of  the  core  storage  required  by  NETSIDE. 
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Problem  Number 


Number  of  Nodes 


Number  of  Arcs 


Transportation  Problems 

1 

100  X  100 

2 

100  X  100 

3 

100  X  100 

4 

100  X  100 

5 

100  X  100 

6 

150  X  150 

7 

150  X  150 

8 

150  X  150 

9 

150  X  150 

10 

150  X  150 

Assignment  Problems 

11 

200  X  200 

12 

200  X  200 

13 

200  X  200 

14 

200  X  200 

15 

200  X  200 

Capacitated  Network  Problems 

16 

400 

17 

400 

18 

400 

19 

400 

20 

400 

21 

400 

22 

400 

23 

400 

24 

400 

25 

400 

26  .  :  _ 

400 

27 

400 

Uncapacitated  Network  Problems 

Table  2  Comparison  of  NETSIDE  and  EQFLO  on  35  Test  Problems 
(Same  step  size  for  every  equal  flow  pair) 


Problem 

Number 

NETSIDE 

EQFLO 

Optimal 

Objective 

Total 

Time 

Time 

Z  of  Optimal 
at  Termination 

Total 

Pure 

Network 

Lower 

Bound 

Upper 

Bound 

Lower  Upper 
Bound  Bound 

1 

2694547 

30 

50 

7 

19 

24 

98 

108 

2 

2350637 

24 

58 

7 

23 

28 

95 

105 

3 

1939836 

27 

127 

9 

55 

63 

99 

110 

4 

1612265 

33 

79 

10 

33 

36 

98 

108 

5 

1480741 

33 

40 

12 

13 

15 

97 

108 

6 

2472907 

71 

46 

22 

9 

15 

98 

108 

7 

2236784 

96 

59 

28 

14 

17 

97 

107 

8 

2223900 

84 

58 

32 

11 

15 

99 

108 

9 

1839835 

115 

44 

36 

6 

2 

98 

105 

10 

2291942 

105 

79 

36 

19 

24 

96 

106 

11* 

4992 

135 

51 

17 

11 

23 

98 

108 

12 

3573 

105 

93 

23 

20 

50 

95 

105 

13 

3142 

103 

78 

27 

16 

35 

98 

108 

14 

2787 

118 

34 

31 

1 

2 

99 

101 

15 

2795 

150 

127 

35 

31 

61 

97 

108 

16 

82161432 

43 

8 

6 

1 

1 

99 

107 

17 

45601025 

66 

13 

8 

4 

1 

99 

105 

18 

81600312 

40 

8 

6 

1 

1 

99 

106 

19 

45601025 

66 

12 

8 

3 

1 

99 

102 

20 

74065202 

40 

9 

6 

2 

1 

99 

108 

21 

40137087 

44 

11 

8 

2 

1 

99 

101 

22 

73429862 

32 

8 

6 

1 

1 

99 

109 

23 

39354594 

33 

11 

3 

2 

1 

99 

101 

24 

85926653 

91 

7 

3 

3 

1 

98 

104 

25 

58203746 

66 

9 

5 

3 

1 

99 

101 

26 

74267081 

65 

6 

3 

2 

1 

97 

102 

27 

47295659 

57 

7 

4 

2 

1 

99 

107 

28 

131316225 

201 

31 

20 

9 

2 

99 

107 

29 

113594497 

260 

167 

25 

72 

70 

98 

107 

30 

90569484 

337 

243 

23 

111 

109 

91 

106 

31 

34943754 

296 

44 

24 

16 

4 

99 

109 

32 

180390305 

529 

80 

48 

25 

7 

98 

109 

33 

205246112 

453 

83 

47 

23 

13 

98 

108 

34 

166247998 

477 

95 

51 

24 

20 

96 

106 

35 

*  Pn 

163964307 

-  10. 

503 

68 

52 

11 

5 

99 

107 

30 


Table  3  Comparison  of  NETSIDE  and  EQFLO  on  35  Test  Problems 
(Different  step  sizes  for  different  pairs) 


Problem 

Number 

NETSIDE 

_  . .  - - 

- 

EQFLO 

Optimal 

Total 

Time 

I  of  Optimal 

Objective 

Time 

at  Termination 

Total 

Pure 

Lower 

Upper 

Lower  Upper 

Network 

Bound 

Bound 

Bound  Bound 

1 

2694547 

30 

47 

7 

16 

24 

97 

108 

2 

2350637 

24 

50 

7 

19 

24 

94 

105 

3 

1939836 

27 

105 

9 

42 

54 

99 

109 

4 

1612265 

33 

74 

10 

29 

35 

98 

108 

5 

1480741 

33 

37 

12 

10 

15 

95 

106 

6 

2472907 

71 

46 

22 

9 

15 

98 

107 

7 

2236784 

96 

57 

28 

13 

16 

97 

105 

8 

2223900 

84 

42 

32 

4 

6 

98 

109 

9 

1839835 

115 

44 

36 

6 

2 

98 

105 

10 

2291942 

105 

76 

36 

19 

21 

96 

105 

11* 

4992 

135 

53 

17 

9 

27 

94 

104 

12 

3573 

105 

79 

23 

14 

4* 

95 

105 

13 

3142 

103 

63 

27 

11 

25 

98 

108 

14 

2787 

113 

34 

31 

1 

2 

99 

101 

15 

2795 

150 

108 

35 

24 

49 

98 

107 

16 

82161432 

43 

8 

6 

1 

1 

99 

107 

17 

45601025 

66 

13 

8 

4 

1 

99 

105 

18 

81600312 

40 

8 

6 

1 

1 

99 

106 

19 

45601025 

66 

12 

8 

3 

1 

99 

102 

20 

74065202 

40 

9 

6. 

2 

1 

99 

108 

21 

40137087 

44 

11 

8 

2 

1 

99 

101 

22 

73429862 

32 

8 

6 

1 

1 

99 

109 

23 

39354594 

33 

11 

8 

2 

1 

99 

101 

24 

85926653 

91 

7 

3 

3 

1 

98 

104 

25 

58203746 

66 

9 

5 

3 

1 

99 

101 

26 

74267081 

65 

6 

3 

2 

1 

97 

102 

27 

47295659 

57 

7 

4 

2 

1 

99 

107 

28 

131316225 

201 

31 

20 

9 

2 

99 

107 

29 

113594497 

.260 

81 

25 

28 

28 

97 

107 

30 

90569484 

337 

148 

23 

62 
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Table  4  Summary  of  Computational  Results 
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